J* 



T\ MASSACHUSETTS INSTITUTE OF TECHNOLOGY 

ARTIFICIAL INTELLIGENCE LABORATORY 

and 

CENTER FOR BIOLOGICAL INFORMATION PROCESSING 
WHITAKER COLLEGE 



A.I. Memo 786 November 1984 

C.B.I.P. Paper 010 



The Coordination of Arm Movements 
An Experimentally Confirmed Mathematical Model 



Tamar Flash and Neville Hogan 



Abstract: This paper presents studies of the coordination of voluntary human arm 
movements. A mathematical model is formulated which is shown to predict both the 
qualitative features and the quantitative details observed experimentally in planar, multi-joint 
f> arm movements. Coordination is modelled mathematically by defining an objective function, 

a measure of performance for any possible movement. The unique trajectory which yields the 
best performance is determined using dynamic optimization theory. In the work presented 
here the objective function is the square of the magnitude of jerk (rate of change of 
acceleration) of the hand integrated over the entire movement. This is equivalent to 
assuming that a major goal of motor coordination is the production of the smoothest 
possible movement of the hand. The theoretical analysis is based solely on the kinematics 
of movement independent of the dynamics of the musculoskeletal system, and is successful 
only when formulated in terms of the motion of the hand in extracorporal space. The 
implications with respect to movement organization are discussed. 
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/^ Introduction 



How is movement control organized? Which variable(s) are controlled? These questions 
have become a growing concern of motor neurophysiologists [Granit 1981, Stein 1982]. 
Although investigations have traditionally focused on single muscle contractions or single- 
joint movements, these systems cannot reveal the problems confronted by the central 
nervous system in the control of normal multi-joint movements. Even a two-joint motion is 
vastly more complicated than a single-joint motion; in moving from one point to another, 
on what basis does the central nervous system select one specific trajectory from the 
infinite number possible? In what coordinate frame is the trajectory planned? However, 
these complexities offer new research opportunities; investigations of multi-joint movements 
may provide considerable insight into the strategies employed by the central nervous system 
in the control of skilled activities. 

Recently a few studies of the kinematic and dynamic aspects of multi-joint human and 
monkey arm movements have been conducted. The objective of these studies was to identify 
common kinematic features or stereotyped patterns of muscle activation characterizing these 
movements [Morasso 1981, Abend et al. 1982, Soechting & Lacquaniti 1981, Georgopoulos 
rf—s, et al. 1981 , Hollerbach and Flash 1982]. The planning and control of the kinematic aspects of 

arm movements is termed trajectory formation. The term trajectory refers to the configuration 
of the arm in space and to the speed of movement as the hand moves from its initial to its 
final position. Some investigators [Greene 1972, Soechting and Lacquaniti 1981, Saltzman 
1979] have argued that trajectories are planned in joint variables. It has been claimed that 
the CNS uses a strategy of maintaining constant ratios between angular velocities of the 
joints in order to bring about a reduction in the complexity of the control problem by 
reducing the number of degrees of freedom. In contrast to this view, other investigators 
have argued that simplicity of motor control is achieved by planning hand trajectories 
in extracorporal space; joint rotations are then tailored to produce these desired hand 
movements [Lashley 1951, Bernstein 1967]. Recently, this view has gained support from 
studies of planar, unconstrained human and monkey movements [Morasso 1981, Abend et 
al. 1982, Georgopoulos 1981]. When moving the hand between pairs of targets subjects 
tended to generate roughly straight hand trajectories with single-peaked, bell-shaped speed 
profiles; this behavior was independent of the part of the work-space in which the movement 
was performed. Because the common invariant features of these movements were only 
evident in the extracorporal coordinates of the hand, these results are a strong indication 
that planning takes place in terms of hand trajectories rather than joint rotations. 

- i# *" s Can this conclusion be generalised to more complex movements? When subjects were 

instructed to generate curved movements, the single-peaked hand speed profile was not 



^**\ preserved. Although the hand paths appeared smooth, their curvature was not uniform, and 

the hand speed displayed one or more maxima. The minima between two adjacent speed 
peaks corresponded temporally to peaks in the curvature [Abend et al. 1982]. A similar 
temporal coupling between speed and curvature has been observed in handwriting, drawing 
[Viviani & Terzuolo 1980] and in infant reaching movements [von Hofsten 1979]. In this 
paper we will show that this behaviour can be derived from a single organising principle. 
This principle is based on the kinematics of the motion of the hand in extracorporal space, 
and provides further evidence that movements are planned in terms of hand trajectories 
rather than joint rotations. 

To describe this behaviour, a mathematical model of the organization of voluntary arm 
movements is presented. There are, of course, many ways of formulating a mathematical 
description of any given phenomenon. In the work presented here, we have used dynamic 
optimization theory as it permits us to describe an assumed goal of this class of movements 
in a relatively simple formula, and derive from the formula a detailed prediction of the 
kinematics of a large number of specific movements. It will be shown that this mathematical 
model succeeds in accounting for the majority of the kinematic features of planar horizontal 
arm movements described in previous studies [Morasso 1981, Abend et al. 1982]. A number 
of new features of planar horizontal arm movements are also predicted by the model, 
/""Y and some new experimental results confirming these predictions are presented [Flash and 

Hogan 1982]. 

The Mathematical Model 



Briefly, dynamic optimization requires the definition of a criterion function which 
describes the objective of the movement. Generally, this function is expressed mathematically 
as a time integral of a performance index, an algebraic function which may in general 
depend on the system inputs, outputs, and internal variables. A set of differential equations 
are formulated, which describe the response of the system to its inputs. The methods 
of variational calculus and optimal control theory are applied to find the trajectory which 
minimizes this criterion function subject to dynamic constraints imposed by the system 
differential equations and the algebraic constraints imposed at the end points, or during the 
motion. The use of optimization techniques to model natural behavior is appealing because 
of the analogy it bears to the optimization presumed to occur as a result of natural selection. 
Studies of two-joint arm movements have shown that as a result of practice the variability 
in hand trajectories is reduced exponentially with time [Georgopoulos et al. 1981]. The fact 
that only a few executed trajectories emerge, may indicate an underlying adaptive process 
tending to produce movements which optimize certain kinematic or dynamic variables 
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[Abend et al. 1982]. The critical step in the analysis is the selection of an appropriate 
objective function. Here the experimental results are suggestive. Since with learning and 
practice movements tend to be performed more smoothly and gracefully, this may indicate 
an underlying objective of achieving the smoothest movement which carries the hand from 
one equilibrium position to another. This was the point of departure of an optimization-based 
mathematical description of voluntary movements in monkeys [Hogan 1982, 1984, Bizzi et 
al. 1984]. The major qualitative and quantitative features of single-joint forearm movements 
have been successfuly predicted assuming that maximizing smothness may be equated to 
minimizing the mean-square jerk. Jerk is mathematically defined as the rate of change of 
acceleration. The work presented in this paper is the generalisation of that analysis to the 
case of multi-joint motion. 

An important feature of the multi-joint case is that the predicted behavior depends 
critically. on the choice of coordinates in which the criterion function is formulated. Our 
choice of coordinates was again guided by experimental observations, in particular the fact 
that the invariant features of upper limb movements are only evident when hand motion is 
expressed in extracorporal coordinates. 

The position vector of the hand was defined with respect to a laboratory-fixed Cartesian 
coordinate system. Differentiating this position vector three times, Cartesian jerk for the 
hand can be defined. In moving from an initial to a final position in a given time t f , the 
objective function to be minimized is the time integral of the square of the magnitude of 
jerk: 
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x and y are the time-varying hand position coordinates. Mathematical expressions for x(t) 
and y(t) are to be found, which bring the criterion function in equation 1 to a minimum. 

This optimization procedure was used for the description of several experimentally 
observed types of human planar two-joint arm movements: unconstrained point-to- 
point movements, unconstrained curved movements and obstacle-avoidance movements. 
Comparisons of the mathematically predicted trajectories with experimental movement 
records were used to evaluate the success of the model. 



Unconstrained Point-to-Point Movements 

For unconstrained point-to-point movements the objective can be stated as follows: 



^\ "Generate the smoothest motion to bring the hand from the initial position to the 

final position in a given time." 

One might expect the physical system which generates the movement (i.e. the neuro- 
musculoskeletal system) to impose certain constraints on the kinematic or dynamic variables. 
Such constraints might set limits on the hand speed or acceleration, due, for example, 
to limitations on the maximum torques that the system can generate or on how rapidly 
they can be changed [Nelson 1983]. However, none of the movements studied here was 
extremely fast nor did any of the movements require the generation of large forces. As 
the subjects did not operate anywhere near the limits of neuromuscular performance these 
constraints are inoperative. The optimization procedure is presented in Appendix A and 
results in a minimum jerk trajectory which is a fifth order polynomial in time for both x(t) 
and y(t). 

If the constraints imposed by the physical system become important, the problem may 
be solved using the method of Pontryagin [Bryson & Ho 1975]. This method was also 
applied to our problem (Appendix B). The two methods yielded the same expressions for 
the hand trajectory. 

The criterion function determines the form of the movement trajectory. The details 
f*S are determined by the boundary conditions at the onset and termination of the movement. 

Given this information and the duration of the movement, the trajectory of the hand is 
specified in its entirety. No other information is required. Assuming the movement to start 
and end with zero velocity and acceleration, the following expressions for hand trajectory 
are obtained: 

x{t) = x + (x - Xf^lSr* - 6r 5 - 10r 3 ) 

y(t) = yo + (yo - S//)(l5r 4 - 6r 5 - 10r 3 ) l } 

where r = t/t f and x > yo are tn <3 initial hand position coordinates at t = 0. 

Profiles of the predicted hand path (P) and hand tangential velocity (T) are shown in 
Figure 1A. The hand path, P, is illustrated by plotting the predicted hand location every 10 
msec. Ti to T e represent the targets. The plot of the predicted tangential velocity profile 
illustrates the hand speed, V, relative to the maximal speed, V M , as a function of t/t f , 
where t is the time and t f is the movement duration. This trajectory is a straight line between 
initial and final positions with a bell-shaped unimodal velocity profile. Since the predicted 
trajectory depends only on the initial and final positions of the hand, it is invariant under 
translations and rotations. The shape of the predicted trajectory also does not change with 
f m ^-- amplitude or duration of the movement, which merely serve to change the scale of the 

position and time axes respectively. Note also that the expression describing the trajectory 
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Figure 1 (A) A representative example of a predicted minimum-jerk unconstrained point-to-point 
hand trajectory. (B) A representative example of a predicted curved minimum-jerk hand trajectory 
between targets 1 and 4. 



along the line joining the initial and final positions is a single fifth order polynomial in time 
which is identical to the minimum-jerk trajectory for single-joint motion [Hogan 1982, 1984]. 
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/•""n Curved Point-to-Point Movements 

To model curved and obstacle-avoidance movements, it was assumed that the hand 
is required, in its motion between the end points, to pass through a third specified point. 
Hence for these movements the movement objective can be stated as follows: 

"Generate the smoothest motion to bring the hand from the initial position to the 
final position in a given time. The hand must move to the final position through 
a specified point (a via point) at an unspecified time." 

The requirement that the hand should move through a specified via point defines equality 
constraints; i.e., the hand position is prescribed at some intermediate time, t u between the 
onset and the termination of the movement. Hence, if the location of the via point with 
respect to a Cartesian coordinate system is given by the coordinates (xi, yi), the equality 
constraints on the hand position coordinates, x(t) and y(t), at time t lt are: 

x{h) = Xi 
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The time <i at which the hand has to pass through this point is not a priori specified but 
is derived from the optimization procedure. It is one of the predictions of the mathematical 
model. Problems of this kind are known as dynamic optimization problems with interior 
point equality constraints and techniques have been established for their solution [Bryson 
& Ho 1975]. The optimal trajectory for the entire movement between the initial and final 
position is derived in appendix C. For the particular problem treated here the technique 
used requires the continuity of the velocity and acceleration at this intermediate time. Since 
discontinuities of hand velocities and accelerations would require infinite accelerations and 
jerks, this is a reasonable requirement for arm movement. The details of the predicted 
movements, including the time at which the hand passes through the via point, depend 
on the chosen boundary conditions at time t — and t = t f and on the specified interior 
point. 

Applying the optimization technique one obtains an expression for the position 
component x(t) at all times t < t x : 



x-( r ) = iLLJr^lbr* - 30r 3 ) + r?(80r 3 - 30r 4 ) - 60r 3 r? + 3O7S - 6r 5 ) 

+c 1 fl5r 4 -10r 3 -6r 5 ]]+a:o (4) 



f\ and for times t > t x the expression is: 

x + {t) = t4( ^if r i( 15T< ~ 30r3 + 30t - 15 ) + *1(-30r 4 + 80r 3 - 60r 2 + 10) ) 
+ci( -6r 5 + 15r 4 - 10r 3 + 1 ) J + x f 
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These equations depend on the dimensionless variables n = h/tj and r = t/t f . Similar 
expressions are obtained for y + {r) and for y~{r) with constant coefficients, n 2 and c 2 
replacing ni and ci in equations 4 and 5. These four coefficients, wx, c 1( tt 2 and c 2 which 
depend on the position coordinates at the boundaries and at the interior point and on t x 
are defined in Appendix C. 

A typical predicted minimum-jerk curved movement is shown in Figure 1B. Displayed 
are the hand path (P) and the hand tangential velocity (T) and curvature (C) profiles. In 
addition, the two velocity profiles along the x axis, V x , and the y axis, V y , are also shown 
in Figure 1B. The meanings of d it d 2 and d 3 are explained in the text below. As this 
S*\ figure demonstrates, the predicted hand curvature increases, reaches a maximum, and then 

decreases again as the movement approaches the final position. The tangential velocity 
curve exhibits two peaks and the minimum in the tangential velocity curve occurs at the 
same time as the maximum curvature point. Again, the predicted trajectories are invariant 
under translation, rotation, time scaling and amplitude scaling. 



Experimental Methods 



To evaluate the model, simulated trajectories were compared to measured hand 
trajectories. Planar horizontal arm movements were recorded, using the apparatus shown 
in Figure 2. The subject was seated and held the handle of a two-link mechanical 
manipulandum. The shoulder was restrained throughout these experiments. The subject's 
wrist was braced in some of the trials, in others it was free to move. The subject was 
instructed to move the handle of the manipulandum under targets mounted on a horizontal 
plexiglass panel above the apparatus, upon their illumination. Visual information about the 
arm location was eliminated by darkening the room in which recordings took place. The joint 
angles of the mechanical apparatus were monitored by means of precision potentiometers. 
f ^ The potentiometer voltage signals were digitized at a rate of 100 samples per second and 

stored on magnetic tapes. These signals were filtered with an upper cut-off frequency of 5.2 
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Figure 2 Experimental apparatus for measuring arm trajectories in a horizontal plane. Ti to r 6 are 
the LED targets. 61 and 2 are respectively the subject's shoulder and elbow joint angles. Movement 
of the handle was measured by way of the potentiometers, P1 and P2. 

Hz. A calibration procedure, executed before each experiment, determined the parameters 
for the conversion of these voltage signals to apparatus joint angles. Off-line data analysis 
was performed to compute the subject's hand position, and, based on the measured 
geometry of the subject's upper extremity, joint angles were also computed. Lagrange 
polynomial differentiation was used to obtain joint angular velocities, hand velocities and 
the curvature of the hand path. Handle movements as small as 1.0 mm could be detected 
and the error in the computed hand velocity was less than 4%. The error in computed hand 
curvature was less than 6%. 

This apparatus and procedure are the same as used by Abend et al. [1982] and in 
this paper their results are compared to predictions of the mathematical model. The model 
also predicts several additional features of movement and to test these predictions new 
experiments were performed. Abend et al. used four experimental paradigms. In the first 
paradigm the subjects were instructed to move the hand from one target to another upon 
illumination of the second target. No instructions were given regarding the path between 
the two targets. We used the same experimental paradigm instructing the subjects to move 
at various speeds by specifying movement durations ranging from about 0.5 to 1.0 sees. 
In the second paradigm subjects were instructed to move the hand to a target along 
a self-generated curved path. No guide or obstacle was presented. In the third group 
of experiments, subjects were instructed to follow curved guides. These were constant 



j0*\ curvature arcs placed on top of the plexiglass panel. The fourth paradigm required obstacle 

avoidance. Subjects were instructed to move their hand from one target to another while 
avoiding an "obstacle" represented by a row of red light-emitting diodes. The maximum 
length of obstacles used in these experiments was 25 cm. The obstacle was placed between 
the targets with its long axis perpendicular to the line connecting the two targets. 

In the above paradigms the subjects were not constrained to move the hand through 
any specific intermediate point on route to the target. We were therefore interested in 
determining what effect the introduction of a real via point would have on the movement. 
Experiments were conducted by us in which subjects were instructed to generate continuous 
movements from one target to another through an intermediate target. The intermediate 
target was lit throughout the movement. Movements involving the same set of three targets 
were performed at the subject's preferred speed first and then at a faster speed. Overall 
nine subjects participated in these experiments. 

In most experiments the movements were confined to elbow and shoulder rotations 
in a horizontal plane at the subject's shoulder level. In some experiments, however, wrist 
movements were also allowed and/or subjects were instructed to move their arm in a 
horizontal plane passing through their waist. 
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To determine model-predicted trajectories for comparison, the values of several 
kinematic variables were derived from the experimental data. For straight trajectories, 
these values were the initial and final target locations and the movement durations. For 
unconstrained curved and obstacle avoidance movements, the location of the via point also 
had to be specified. This was determined from the location of the maximum curvature point. 
In experiments with an intermediate target, the location of the via point was taken to be 
the location of this target. These parameters were substituted in the polynomials derived 
from the dynamic optimization and hand positions were computed every 10 msec. These 
positions were then differentiated to derive hand velocities and differentiated again to derive 
hand accelerations and the curvature 1 was computed. 



'Hand speed (T) and trajectory curvature (C) are defined as: T = \/(i) 2 + (y) 2 , and C = *$ ■ ?* . 
where i, y are the hand velocity components and '£, y' are the acceleration components. 
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Figure 3 Overlapped predicted (solid lines) and measured (dotted lines) hand paths (a), speeds 
(b) and acceleration components along y axis (c), and along x axis (d) for two unconstrained 
point-to-point movements. (A) A movement between targets 3 and 6. (B) A movement between 
targets 1 and 4. 

Results. 
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Typical experimental results for two unconstrained point-to-point movements are 
shown in Figure 3 superimposed on the predicted minimum-jerk trajectories for the same 
movements. As this figure shows, there is a good qualitative and quantitative match between 
the predicted and actual trajectories. There is a good agreement between the predicted and 
measured steepness of the rising and falling parts of the tangential velocity and acceleration 
curves, and the time at which maximum acceleration is reached. The predicted trajectory 
does not capture the asymmetry in the tangential velocity profile of the measured trajectory. 
The difference, however is quite small. There are also some slight discrepancies between 
the measured and predicted acceleration profiles. 

Typical examples of the fit between the temporal behavior of real and simulated 
shoulder and elbow angles and angular velocities are shown in Figure 4. To obtain joint 
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Figure 4 Measured and predicted shoulder, S, and elbow, E, joint angles and angular velocity 
profiles. (A) A movement between targets 3 and 6. (B) A movement between targets 2 and 5. 

angles and angular velocities, the inverse kinematics problem [Brady et al, 1982] was solved. 
Since for horizontal planar two-joint movements hand trajectories uniquely define angular 
trajectories, this computation serves only to permit an additional comparison between the 
theoretical and experimental results. As can be seen in Figure 4 the actual and predicted 
angular velocity profiles agree quite well and there exist only minor discrepancies between 
the simulated and real trajectories. 

A typical example for a predicted curved trajectory was shown above in Figure 1B. As 
this example demonstrates, the predicted trajectories display the majority of the observed 
features of real curved hand trajectories as discussed by Abend et al. [1982]. If real curved 
hand trajectories can be described adequately by the proposed model, they should display 
the kinematic features exhibited by the simulated trajectories. To facilitate the explanation 
of the nature of these predictions a simple geometric construction, shown in Figure 1B, is 
used. The movement end-points are connected by a straight line. The perpendicular to this 
line which passes through the via point is drawn! This second line, of length d 3 , divides 
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Figure 5 Representative examples of comparisons between measured (Real) and predicted (Model) 
trajectories from a "via-point" experiment in which the movement reversed its direction along the 
y axis. (A) The intermediate target was located at equal distances from the initial and finai targets. 
(B) The intermediate target was at equal distances from the movement end-points but closer to the 
line connecting them. 

the line connecting the two movement end-points into two segments of lenghts di and d 2 . 
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The first predicted kinematic feature relates to the dependency of the shape of the 
hand velocity and curvature profiles on the location of the via point. If the via point is 
displaced in any direction, towards either the initial or the final target, so that one of the 
two segments d x or d 2 is longer, the amplitude of the velocity peak on the corresponding 
portion of the movement will be higher. The second predicted kinematic feature relates to 
the depth of velocity valley, and the height of the curvature peak. For more highly curved 
movements the tangential velocity dip is more pronounced. Hence, for two movements of 
the same duration with the same values of d x and d 2 , the amplitude of the hand curvature 
peak will be higher and the tangential velocity valley will be deeper for the movement with 
the larger value of d 3 . 
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Figure 6 Representative examples of comparisons between measured (Real, left columns) and 
predicted (Model, right columns) trajectories from a "via-point" experiment in which the movement 
reversed its direction along the y axis, (A) The intermediate target was closer to the initial target. 
(B) The intermediate target was closer to the final target. 

To test the validity of these two predictions, quantitatively simulated and measured 
hand trajectories were compared for movements through intermediate targets, unconstrained 
curved movements and obstacle-avoidance movements. In all the following figures measured 
hand paths (P), hand speed profiles (T), hand curvature profiles (C) and the profiles of the 
two velocity components along the two orthogonal axes (V x and V y ) are displayed in the left 
column (Real) and the corresponding plots for the minimum jerk trajectories are displayed 
in the right column (Model). Comparisons of predicted and measured hand trajectories for 
arm movements generated by subjects in experiments with intermediate targets are shown 
in Figures 5-7. 
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Figures 5 and 6 describe results from trials in which the location of the via point, 
with respect to the initial and final targets, required the reversal of the direction of the 
hand movement along the y axis. In this group of movements, the intermediate target 
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Figure 7 Representative examples of comparisons between measured (Real) and predicted (Model) 
trajectories from a "via-point" experiment in which the movement reversed its direction along both 
x and y axes. 

was first placed at equal distances from the initial and final targets, closer to (Figure 5A) 
or further away (Figure 5B) from the line connecting these two targets, and then was 
symmetrically offset from the center in both directions along a line passing parallel to the 
y axis (Figure 6A and Figure 6B). As these figures show there is a good qualitative and 
quantitative fit between the measured and predicted trajectories. Furthermore, the measured 
hand trajectories validated the above predictions. 

Figure 7 shows results from trials in which the locations of the intermediate targets 
required a reversal of movement direction along both the x and y axes. As this figure 
demonstrates, the fit of the predicted to the real trajectories was equally good under 
conditions of both translation and rotation of the movements. 

Two obstacle-avoidance movements are compared to the corresponding simulated 
trajectories in Figure 8. Two unconstrained curved movements are compared to simulated 
trajectories in Figure 9. Again all of the above predicted features are exhibited. Note also the 
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Figure 8 Representative examples of comparisons between measured and predicted trajectories 
for two obstacle-avoidance movements. 



similarities in the kinematic characteristics of obstacle avoidance movements, unconstrained 
movements, and movements through intermediate targets. 

Additional predictions can be made from the model. These are the "isochrony principle" 
and the scaling of the hand trajectories with speed. The isochrony principle states that the 
times it takes the hand to move along the two portions of the movement, from the start 
to the via point, and from the via point to the end of the movement, should be roughly 
equal, except for cases in which the via point is very close to either one of the movement 
end-points. Hence if ti/(ti + t 2 ) values are plotted versus di/(d 1 + d 2 ) values, for several 
different movements, *i/(*i + t 2 ) should be expected to be roughly equal to 0.5 for most 
of di/(di + d 2 ) values between and 1 except for d x j{d x + d 2 ) values quite close to or 
1. It is also predicted that for similar values of rfi/(di + d 2 ), the values of <i/(ii + t 2 ) will 
be independent of d 3 . This behavior of the minimum-jerk trajectories is consistent with 
the differences in heights of the two velocity peaks on either sides of the velocity valley 
because, if the hand travels along both movement segments in roughly the same time, then 
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Figure 9 Representative examples of comparisons between measured and predicted trajectories 
for two unconstrained curved movements. 
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the velocity amplitude for the longer segment should be higher. 

In Figure 10A, points corresponding to measured ti/(ti + t 2 ) values for various 
experimentally set d 1 /(d 1 + d 2 ) parameters, are superimposed on the predicted *i/(ti +^2) 
versus c?i/(rfi + d 2 ) curve. The measured *i/{*i + t 2 ) values were derived from both 
unconstrained curved movements and movements through intermediate targets. Similarly, in 
Figure 10B points corresponding to measured t 1 /(t 1 + t 2 ) versus d s for various movements 
with similar values of di/(di +d 2 ), are superimposed on the corresponding predicted curve. 
In both figures the real points fall close to the predicted curves; within experimental errors 
the measured trajectories behave as predicted by the isochrony principle. 

The predicted curved trajectories, like the straight trajectories, scale with time as 
the movement duration changes. Therefore, the path of the hand should be the same 
for fast and slow movements, while the hand speed throughout the movement should 
scale numerically by the ratio between the fast and slow movement durations. The scaling 
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Figure 10 (A) A plot of the predicted ti/(h +t 2 ) values versus di/(di + d 2 ) (solid line). Superimposed 
are measured ti/{ti+t 2 ) values derived from unconstrained curved motions (circles), and from 
motions through intermediate targets (squares). (B) A plot of the predicted ti/(ti + t 2 ) values versus 
d 3 (solid line). Superimposed are measured values derived from unconstrained curved motions 
(circles), and from motions through intermediate targets (squares). 
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prediction is validated by the measured movements as demonstrated in Figure 11 A, in which 
the movement duration is 1.1 sec versus Figure 11B, with a movement duration of 0.85 sec. 
Moreover, as can be seen from this figure, the fact that the subjects move faster or slower 
does not affect the model performance. 
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Figure 11 A representative example of the scaling with speed of curved trajectories. (A) A 1.1 
sec. movement via an intermediate target. (B) A 0.85 sec. movement through the same targets. 

The model predicts the same hand trajectories regardless of the specific joints involved 
in the generation of the movement or whether they are generated in a horizontal plane at 
the level of the subject shoulder or at other levels. This was found to be the case in the 
work of Abend et al. [1982]. Similarly, in our own experiments, in trials in which the hand 
brace was removed and rotations at the wrist joint were allowed, and in trials in which the 
movements were generated at a horizontal plane passing at the level of the subject 's waist, 
the trajectories were no different from the trajectories described in the above figures. 

The model succeeds in accounting for the majority of the kinematic features of the real 
trajectories. A few fine-grained details of the real movements, however, were not captured 
in the simulated trajectories. These include, for example, the tendency of the first velocity 
peak to be higher than predicted, in movements between targets 1 and 4. Another example 
is the existence of irregularities in the hand path (the little "hooks") and hand speed as 
the hand approaches the goal target (e.g. Figures 7, 8). 
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^ Discussion. 



The mathematical model presented in this paper matches observed human planar 
two-joint arm movements. The organisation of the movements was modelled through the 
criterion function. From the optimization of the criterion function explicit analytic expressions 
for the description of many different hand trajectories were derived. The derivation of these 
expressions depended on the specification of a small number of parameters; in the case 
of unconstrained point-to-point movements, only the movement duration and the position, 
velocity and acceleration at the end-points were required. Application of the same model 
to curved movements required only the addition of the position coordinates of a via point. 
In return, the model yielded the time at which the hand passes through the via point and 
a detailed time-history of the hand positions, velocities etc. The mathematical analysis, 
particularly the dynamic optimization methods used, bring about a dramatic reduction in 
the dimensionality of the problem of describing movement kinematics. 

The mathematical model presented in this paper cannot be dismissed as merely 
elaborate curve fitting as it leads to several testable predictions: invariance of hand 
trajectories under translation and rotation in the workspace, scaling of hand trajectories with 
a0**^ time and amplitude, temporal coupling between curvature and speed, and the "isochrony 

principle". All of these kinematic characteristics have been observed in the experiments 
reported in this paper and in the motor control literature. For example, it has repeatedly been 
observed that peak velocity increases as an approximately linear function of the distance 
to be traveled, so as to keep the duration of the movement roughly constant [Viviani & 
Terzuolo 1982]. Temporal coupling between hand curvature and speed has been previously 
described [Viviani & Terzuolo 1980, Abend et al. 1982] and movement scaling with time 
has also been observed for many different movements including two-joint arm movements, 
three dimensional reaching movements, handwriting, wrist movements etc. [Schmidt 1980]. 

How does the mathematical model presented here compare with alternative descriptions? 
Only a few models for the mechanisms underlying trajectory formation have been proposed. 
Given the tendency to generate roughly straight hand paths in point-to-point movements 
and the phenomenon of the temporal coupling of hand curvature and speed, Abend et 
al. [1982] and Viviani and Terzuolo [1982] have suggested that this phenomenon is due 
to a central mechanism which plans trajectories as sequences of movement segments 
which are then overlapped in time. Viviani and Terzuolo [1982] claim that movements obey 
an "isogony principle"; the ratio between tangential velocity and radius of curvature is 
ptecewise constant. Our analysis of arm movements did not support this finding. Morasso 
^^ and Mussa-lvaldi [1982] suggest that curved movements are generated from separate 

strokes, each stroke characterized by various geometric parameters, such as length, total 
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f\ angular change etc. Different movements may then be composed from different strokes, 

with different geometric parameters. Their model implies a need for storage and retrieval 
of strokes from a dictionary in memory; This dictionary would have to be very large. 

These descriptions, however, do not make clear what rules and principles govern 
the choice of kinematic variables, nor why there is a tendency to generate straight paths 
in unconstrained movements in the first place. For example, the work by Morasso and 
Mussa-lvaldi [1982], which suggests a computational model for trajectory planning, does 
not clarify what algorithm is used to overlap succeeding strokes in time nor does it suggest 
any explanation for the isochrony principle. By contrast, based on dynamic optimization of a 
single criterion function the model presented here accounted for the entire hand trajectory 
in both unconstrained and curved point-to-point movements. Furthermore, this model could 
be extended to describe more complicated movements by using more via points. However, 
a segmentation mechanism in the execution of hand trajectories is not ruled out by the 
optimization-based model and it is possible that such a mechanism underlies the planning 
of long and complicated sequences of movement. Whatever the planning process, the work 
presented here makes it clear that not every inflection point in the hand path nor every 
velocity valley necessarily imply the beginning of a new stroke. The success of the model 
presented here indicates that a more basic principle than mere segmentation underlies the 
/ P ** N specification of hand trajectories. Some evidence in support of this suggestion can be found 

in Wing's [1978] studies of handwriting in which it was shown that there exists a positive 
temporal correlation between the up and down strokes in the letter n, for example. This 
may indicate that an up-down stroke is planned as a single unit rather than two separate 
strokes. 

This work also offers a plausible algorithm for Cartesian end-effector trajectory planning 
for artificial manipulators. Current methods of trajectory planning in robotics usually constrain 
the end-effector to follow straight lines using first order polynomials [Taylor 1979]. For 
more complicated movements or movements through via points, end-effector trajectories are 
planned as series of simple second-order polynomials which are smoothly joined together 
[Brady 1982] at knot points. While this method explicitly plans end-effector trajectories and 
assures that the end point will stay within the workspace it does not guarantee low values 
of jerk. By contrast the jerk minimization procedure presented here provides a rigorous 
way of planning trajectories of motions between equilibrium points and motions through via 
points while guaranteeing predictable and well-behaved trajectories. 

The results presented here are independent of the physical system which generates 

the motion. No representation of the arm dynamics is required, nor is it necessary to make 

f"*^, any assumption about the form of the neural input signals as in some other proposed 

models which attempt to account for movement kinematics [Zangenmeister et al. 1981]. 
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/"■n Hence, the success of the model presented here supports the theories which view the 

motor system as being divided between higher levels which plan ideal trajectories for the 
end-effector, and lower level processes which translate them into torques and forces. These 
theories suggest that at higher levels there exists a kinematic representation of movement 
which does not take into account the mechanical nature of the actual effectors [Bernstein 
1967]. 

The optimization model offers a new insight into these theories: the minimum-jerk 
movement is independent of the neuromuscular dynamics only if the demands of the 
movement lie within the capabilities of the neromuscular system. If one (or more) of the 
neuromuscular performance limits are reached, they impose a constraint On the achievable 
movements and the planning and execution processes cannot be clearly separated. However, 
dynamic optimization theory is sufficiently general to cover the case of performance-limited 
movements and may permit a prediction of the interaction between movement kinematics 
and neuromuscular dynamics as the limits of performance are approached, for example 
during fast movements. 

During the movement of a multi-joint limb, the generation of appropriate joint torques 
for a trajectory is complicated by the presence of significant joint interactions due to 
inertia!, centripetal, and.coriolis forces [Hollerbach and Flash, 1982]. Since these interaction 
effects are not present in the case of single-joint movements, strategies for the solution of 
dynamic problems for single-joint movements may not generalize to multi-joint movements. 
By contrast, the minimization of mean-squared jerk has been applied successfully to both 
single-joint [Hogan 1982, Hogan 1984] and multi-joint movements [Flash 1983]. By an 
appropriate choice of the boundary conditions (i.e the acceleration at the end of the motion 
is not constrained to be zero) the minimim-jerk model can also produce an excellent 
fit to observed repetetive movements [Nelson 1983]. Similarily, curved and straight hand 
movements were predicted here from the same criterion function while at the dynamic level 
each different movement requires very different patterns of joint torques. This indicates that 
a single unifying principle underlies the planning and the kinematic representation of all 
these movements. These results also indicate that at the higher levels a small number of 
general principles may be applied to the planning and coordination of movement kinematics, 
whereas at the lower levels less general strategies may be used to compensate for the 
difference in the dynamics of multi-jointed and single-jointed systems. 

It is important to note that in the mathematical analysis presented here the minimum-jerk 
criterion function is expressed as the sum of the squares of the third derivatives of the 
coordinates of the hand in the extracorporal space. If the criterion function were expressed 
in joint coordinates (as the sum of squares of the third time-derivatives of the joint angles) the 
observed kinematic features and invariances in the hand trajectory would not be predicted. 
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f*\ Consequently, this work provides strong support for the hypothesis that movements are 

planned in terms of the motion of a "disembodied" hand moving in extracorporal space. 

The concept of hierarchical motor planning [Saltzman 1979, Keele 1981] suggests 
also that same general and abstract internal representation of movement is used each time 
a movement is generated with temporal and spatial parameters chosen for that specific 
movement. Movement durations and spatial position cues were suggested as such specifiers 
[Keele 1981]. This would allow the CNS to use the same reference spatial coordinates for 
coding both visual information and motor actions and for learning, storage and retrieval 
of informaion in other skills such as drawing or handwriting. The minimum-jerk model is 
completely consistent with this hypothesis. The information necessary to obtain a detailed 
specification of the trajectory of the hand between the two targets was the movement 
duration and the extracorporal (e.g. Cartesian) locations of the targets. 

The" fact that hand trajectories are seldom performed twice in exactly the same way 
can be used to argue that different rules are used each time a movement is generated. 
Such variability in hand trajectories generated by different subjects or the same subject 
on succeeding trials may however be due to slight changes in the perceived locations of 
points in space to which the hand moves, or through which it passes. Another possible 
source of movement variability might be motor variability, since the final product depends 
on the translation of the motor plan into muscular activity and on the interaction of the arm 
with the environment [Keele 1981]. 
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Accordingly, one physiological interpretation for the qualitative and quantitative fit 
of the experimentally recorded movements to minimum-jerk movements is that the CNS 
explicitly uses a trajectory planning strategy, which is captured by the mathematical model 
presented here, in order to translate task objectives, coded by exteroceptive position cues, 
into trajectory plans. An alternative explanation is that the smoothness of hand trajectories 
is an outcome of the intrinsic properties of the neural and muskuloskeletal hardware. This 
alternative explains the observations at the level of the neuromuscular structure rather 
than at the level of higher cognitive processes, but the basic principle is the same: that 
evolutionary adaptations have led to an optimization of biological movements. However, 
as mentioned above, the model presented here is based on the jerk of the hand in 
extracorporal coordinates and is unsuccessful when expressed in joint angles. Therefore, 
it can be argued that the smoothness of contraction of individual muscles cannot be the 
cause of the agreement between theory and experiment. 

The rationale for jerk minimization in biological trajectory planning does not lend itself 

to self-evident, casual explanations. Given the fact that the movements under consideration 

f~^ occur at moderate speeds and do not subject the system to undue stress, it is unlikely that 

such a strategy has evolved to minimize the "wear and tear" on the system. It is possible 
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that the objective is to minimize unwanted, abrupt changes in the forces transmitted to 
objects carried by the hand. Another possibility [Hogan 1984] is that the objective is to 
maximize the predictability of the trajectory, which is consistent with minimizing its higher 
time-derivatives. To discriminate between these and other possiblities will require further 
work. 

Finally, it is not suggested that minimization of jerk is the single objective underlying 
all movements. Minimization of mean-squared jerk is a mathematical model of one 
movement objective, the production of smooth, graceful movements. Alternatives such as 
the minimization of mean-squared acceleration and mean-squared snap (the fourth time- 
derivative of position) have been explored [Flash 1983]. This analysis has demonstrated that 
the minimization of mean-squared acceleration results in parabolic tangential velocity profiles 
and non-zero accelerations at the movement end-points, unlike the actual movements. On 
the other hand minimization of mean-squared snap provides a good fit to the observed 
experimental data. The limited resolution of the experimental data, however, has not allowed 
us to establish unequivocally which one of the two models, jerk minimization or snap 
minimization, offers a better fit. Further pursuit of the question of what objectives are 
optimized in human movements may help in clarifying what other principles underly motor 
planning [Pew and Baron 1978]. 
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The optimization procedures used in this work are described in the following appendices. 
A comprehensive introduction to dynamic optimization and variational calculus is found in 
Bryson and Ho [1975] and Pontryagin et al. [1962]. 



Appendix A: Unconstrained Optimization. 

Unconstrained point-to-point movements. 

We want to minimize the following objective function: 

Generally, for any function x(t), which is sufficiently differentiable in the interval <t <t f , 

and for any performance index L[t, x,x,'x, *£f] which is integrable over the same interval, 

the unconstrained cost function: 

C(x{t))=J L[t,x,x,x — ]dt {Al) 

assumes an extremum when x(t) is the solution of Euler-Potsson equation: 

at dt ( dx } + [ l) dt»d(»x)-° [A2) 

Since in our case 

d 3 .dx*. d 3 ,af. n ltA . 

We can uncouple the terms depending on the two position components to get: 

The resulting solution to differential equations of this kind is given by a fifth order polynomial: 

x(t) — clq + a x t + a^t 2 + a^t 3 4- a^t* + ast 5 



we get the equation 



y (t) = b + bit + b 2 t 2 + b 3 t 3 + b^ + b 5 t 



W) 
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j/r*\ Appendix B: Dynamic Optimization. 

The optimization method. 

Generally, optimization problems similar to the problem solved here, involve a system 
which can be described by a set of nonlinear differential equations: 

k = f{s(t),u{t),t} (Bl) 

where s{t) is a n vector function of state variables and u(t) is a m vector control 
function. The problem is to find the control u(t) which in carrying the system from an initial 
state «(u) to a final state «(</)> the cost function C(t) is optimized. C(t) is defined as: 

C[t)~ J L[s{t),u(t),t]dt (52) 

Jo 

where L[s{t), u(t), t] is the performance index. 

This problem can be solved by the method of Pontryagin [Pontryagin et al. 1962]. One 
defines a n component co-state (Lagrange multipliers) vector \(t) and a scalar Hamiltonian: 

/"> H[s(t),u{t),t}=L[s{t),u{t),t] + \ T (t)f[s{t),u{t),t} (B3) 

The following differential equations define the necessary conditionns for a minimum to exist: 

m = Mt),u{t),t) {B4) 
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Unconstrained point-to-point movements. 

For our problem we define a state vector z T (t) — [x, y, u, v, z, w] and a control vector 

The components of these vectors are defined by the system equations: 

i = u 
y = v 

u = 'x — z 

• .- (57) 

v = y = w 

z ='&' = jerk z = 5 
u) = y = jcrk y = 7 
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i<*N and the Hamiltonian is: 

H =\ x u + \ y v + X u z + \ v w + \ z 6 + \ w f 



lO 



+ fa* + O (B8) 



The necessary conditions for a minimum to exist are: 



x a 



dt _AV 

dXu, 

The necessary conditions on the control variables are: 



(B9) 



Applying the following boundary conditions 



— = * + x, = o 

— . = 7 + X„, = 


conditions: 




x(0) = x 


2/(0) = J/o 


x{t f ) = X, 


»('/) = V! 


u(0) = 


u{t f ) — 


i>(0) = 


v(t f ) = 


2(0) = 


*(*,) = 


w(0) = 


«)(</) = 



(BIO) 



(Bll) 



We solve these equations and get the following fifth order polynomials for the position 
coordinates x(t) and y(t): 

x{t) = a + a t t + a 2 t 2 + a 3 t 3 + a 4 t 4 + a 5 t 5 
y(t) = b + b t t + b 2 t 2 + b 3 t* + b t t 4 + b 5 t 5 



Appendix C: Via Point Constrained Optimization 

Problems with interior point equality constraints. 

r^ For curved movements we assumed that the hand is constrained to pass through the 

via point at time t x and the coordinates of this point are xj and y v Problems of this kind 
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/•-% are defined as optimal control problems with interior point equality constraints on the state 

variables [Bryson & Ho, 1975], 

For such problems one has a set of constraints at some time t t : 



&(s(h),h) = (CI) 

where K is a q-component vector function. These interior point constraints can be 

augmented to the cost function by a Lagrange multiplier vector a so that the new cost 

function is: 

C = 7L T N+ ['(H-\ T a) (C2) 

Jo 

The solution is obtained by allowing discontinuities in the co-state variables (Lagrange 
coefficients) X(*)'s and in the Hamiltonian H[t,\(t),s(t)}. One can define a vector of Lagrange 
coefficients X + (i) and Hamiltonian H + (t) for t > t x and a vector \~{t) and Hamiltonian H~ 
for t < ti . At time ti these variables satisfy the equations: 

ir(h) = k + (h) + K T : ^- ) (C3) 

^^ . H-M^H+M-z*^ ((74) 



The q components of s. are determined by the constraint equations (C1) while time h is 
fully determined by equation (C4). 

Minimum-jerk trajectories with via point constraints. 

For our specific problem the only constraints are on hand position along both axes at 
time ti; 

y(ti) = vi 
The Hamiltonian H~ for all times, such that t < tx is: 

H~ =X~u~ + \yV~ + \~z~ + \~w~ + \~6~ + X~7~ 

+ ^((7-) 2 +(n 2 ) (<?•) 

and the Hamiltonian H + for times t > ti is: 



f m * s \ 



H + =X+u + + X> + + X+z+ + X+«i+ + X+S+ + Xj 7 + 

+ l((l + ) 2 + (6 + f) (CT) 
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■n*-^ Since the only constraint equations relate to position the only discontinuities are in \ x and 

Xj, and, therefore, according to equation (C3), we get : 

x- - x++», {C8) 

while all the other Lagrange coefficients are continuous at t = t x : 

K(h) = x-(ti) c 

\l(h) = K(h) \' 

Xi(*i) = X-(tx) 

Since time t L is not explicitly specified, the Hamiltonian must be continuous at h as emerges 
from equation (C4) : 

H + (h) = H-fa) (CIO) 



Next we write the necessary conditions for the existence of a minimum. These equations 
are derived separately for t > tx and t < i x . In addition we require continuity of velocities 
and accelerations at t it so that: 

« + (<i) = m~(*i) 

■* (*l) = z (*l) 

These come as a result of the requirements of the optimization technique for continuity of 
the state variables at the interior points [Bryson & Ho, 1975]. 

Solving the equations, defining the necessary conditions for a minimum and applying 
the boundary conditions at t — t and at t — t f , one obtains an expression for x(t) at all 
times t < tx : 

x -(r) = ^( 7r i( r t( 1 ^ 4 - 30r 3 ) + r 3 (80r 3 - 30r 4 ) - &0t 3 t\ + 30^ - 6tA 

+cj 15r 4 - 10r 3 - 6r 5 J J + s (C12) 

and for times t > i x the expression is: 

x+( T ) = -L( Vl f t \(15t* - 30r 3 + 30r - 15) + r 3 (-30r 4 + 80r 3 - 60r 2 + 10) J 
+c t ( -6r 5 + 15r 4 - 10r 3 + 1 J j + a:/ 

-rW + „*&££ «ni) 
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where ci and 7r x are constants, r = t/t/ and 7-1 = *i/t/. 

The same expressions are obtained for y + and for y~ with 7r 2 and c 2 replacing 7r x and 
c! . Substituting :r + (ti) = x~(*i) = xi in equations (C12) and (C13) and solving for ttj and 
a, we obtain the following expressions : 



: ((x f - x )(300r! 5 - 1200n 4 + 1600ti 3 ) 



^rf(l-n) 

+ n 2 (-720i/ + 120ia + 600x ) + (*o ~ *i)(300ti - 200) J (C14) 

Wl = 1 1 ( I/ _ l0 )(l20r 1 5 - 300rj 4 + 200ri 3 ) - 20(x! - x ) J (C15) 

i/ 5 n 5 (l-^) 5V ' 

and for w 2 and c 2 similarly, with y replacing x and y f replacing x f etc . 

Next ,we substitute the expressions for w lt n 2 , c x and c 2 in equation (C10) which 
reduces to 

7ri«(t 1 ) + 7r 2 w(t 1 ) = (C716) 

and we get a polynomial equation in n = t 1 /t f . We find the real roots of this polynomial 
and accept only those roots that lie between and 1. The polynomials we obtained had 
only one acceptable root. We then substitute this value for r x in the expressions for n lt n 2 , 
ci and c 2 and we can finally get the expressions for x(t) and y(t) for the entire movement. 
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